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Abstract 

We consider statistical mechanical properties of the primitive chain network (PCN) model 
for entangled polymers from its dynamic equations. We show that the dynamic equation for 
the segment number of the PCN model does not reduce to the standard Langevin equation 
which satisfies the detailed balance condition. We propose heuristic modifications for the 
PCN dynamic equation for the segment number, to make it reduce to the standard Langevin 
equation. We analyse some equilibrium statistical properties of the modified PCN model, by 
using the effective free energy obtained from the modified PCN dynamic equations. The PCN 
effective free energy can be interpreted as the sum of the ideal Gaussian chain free energy and 
the repulsive interaction energy between slip-links. By using the single chain approximation, 
we calculate several distribution functions of the PCN model. The obtained distribution 
functions are qualitatively different from ones for the simple slip- link model without any direct 
interactions between slip-links. 

1 Introduction 

The primitive chain network (PCN) model [T|[2] is a slip-link type mesoscopic coarse-grained model 
for entangled polymeric systems [3]. To simulate rheological properties efficiently, the PCN model 
represents entangled polymers as a network like structure of which topology dynamically changes. 
Simulations based on the PCN model can reproduce rheological properties of entangled polymers 
well with relatively small computational costs, and it has been applied various systems including 
branched polymers [4, 5 or bidisperse polymers 6 . Although the PCN model achieved success to 
predict various rheological behaviours of entangled polymers, its statistical properties are still not 
fully understood. 

In recent years, various primitive path analysis methods [7M31 have been developed to extract 
statistical properties of network structures in entangled polymers quantitatively from atomistic or 
coarse-grained molecular models (such as the molecular dynamics simulations). It is an interesting 
question whether the extracted statistical data agree with ones obtained by the PCN model, to 
discuss consistency or relation between the models. Quite recently, the network statistics of the 
PCN model were systematically examined [TU]. The statistical properties of the PCN model are 
shown to be qualitatively in agreement with primitive path analysis data. However, the analysis 
is conducted only for simulation results, and from the theoretical view point, how the network 
statistics is determined in the PCN model is still not clear. The PCN model is a dynamical model 
which is constructed in a rather phenomenological way. This makes theoretical analysis for the 
equilibrium statistics of the PCN model difficult. Moreover, the thermodynamic validity of the 
PCN model is not guaranteed from the view point of the statistical mechanics. 

To achieve the thermal equilibrium state, the detailed balance condition is required to be satis- 
fied. If the detailed balance condition is satisfied for a model, we can utilize the standard statistical 
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mechanics to analyze the model. Namely, the equilibrium probability distribution is given by the 
Boltzmann type equilibrium distribution with the thermodynamic potential. In the vicinity of 
the equilibrium state, dynamical behaviours such as the relaxation or linear response properties 
can be related to the correlation functions of fluctuations in equilibrium [TT1113) . The detailed 
balance condition is a strong condition and thus it is not always satisfied in a phenomenologically 
constructed dynamical model. Therefore, it is desired to examine whether the PCN model satisfies 
the detailed balance condition or not. So far, most of the previous works for the PCN model focused 
on dynamical properties such as rheological properties. As far as the authors know, the analysis 
of the PCN model from the view point of statistical mechanics has never been shown explicitly in 
the literatures. (It would be fair to mention that several researchers have noticed that the PCN 
model or some other slip-link models do not satisfy the detailed balance condition, although it has 
not been stressed in published literatures. For some slip-link models, the statistical mechanical 
analysis or the modeling consistent with the detailed balance condition have been done [M l fTS ] .) 

In this work, first we attempt to interpret the dynamic equations of the PCN model as the 
Langevin equations. We show that the PCN dynamic equations do not satisfy the detailed balance 
condition and therefore the PCN model does not have the thermal equilibrium state. Then we 
propose possible modifications to recover the detailed balance condition, in a heuristic way. Even 
if we modify the PCN model, the resulting thermodynamic potential (effective free energy) of the 
PCN model is not identical to the free energy of the slip-link model without any direct interactions. 
In the PCN model, there is the effective repulsive interaction between neighboring slip-links on a 
polymer chain. Based on the obtained effective free energy and the single chain approximation, 
finally we calculate several equilibrium probability distribution functions analytically. We compare 
analytical results with other theories as well as PCN simulation data. 

2 Statistical Mechanical Interpretation of PCN Dynamic 
Equations 

2.1 PCN Dynamic Equations 

In the PCN model [HE], entangled polymers are represented as a network structure which consists 
of nodes (slip-linked points and chain ends) and bonds which connect nodes. For simplicity we 
limit ourselves to monodisperse linear polymer systems with the polymerization index (number of 
segments) being N. The state of the system is described by the set of coarse-grained variables; 
positions of slip-linked nodes and end nodes, numbers of segments between two neighboring nodes, 
and the connectivity information. For convenience, we express the index of the fc-th node in the 
i-th chain as (z, k). We express the position of the (i, k) node as Ri^k, and the number of segments 
between the (z, k—1) and (i, k) nodes as Ni^k- We also express the number of subchains in the z-th 
chain as Zi. The 0th and Z^-th nodes in the i-th chain represent chain ends whereas other nodes 
represent slip-linked nodes. A slip-linked node is assumed to be spatially coupled to the partner 
node. We express this by introducing a connectivity map for the node index, C (the (i, k) node is 
coupled to the C{i, k) node). Two nodes connected by slip-links, (i, k) and C{i, k), share the same 
node position. 

Ri,k = Rc{i,k) (for a-U slip-linked nodes) (1) 

The state of the system can be completely described by {Ri^k}i {Ni^k}, {Zi}, and C. (The 
information of C is not important in the following arguments and thus we do not describe it 
explicitly.) 

In absence of external flow (deformation) fleld, the PCN dynamic equations '2[ are described as 
follows. (Although there are several different versions of the PCN model, in this work we employ 
the version described in Ref[2l) 
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where C is the friction coefficient of a node, ks is the Boltzmann constant, T is temperature, b is 
the segment size, and Ji,k{i) is the flux of the segment number on the (i, k) node. The summation 
in the right hand side of eq ^ is taken for ah nodes connected to the target (topological neighbor 
nodes). The topological neighbor nodes are (j, A: ± 1) and («', k' ± 1) (with («', k') = C(i, k)) for a 
slip-linked node, and (i, fc + 1) or (j, fc — 1) for an end node, pi^k is the segment density along the 
polymer chain on the (i, k) node, and is defined as the arithmetic average of local densities in two 
neighboring bonds. 
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w^^'^ (t) and Wj-^-* (t) are Gaussian white noises which satisfy the following relations. 
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where (...) means the statistical average and 1 is the unit tensor. 

We note that the (osmotic) repulsive force terms [1] are dropped in the dynamic equations ([2])- 
(|4]) . The repulsive force is required only to cancel the artificial attractive interaction in multi chain 
slip-link systems [TB], and its contribution is not essential in the analysis in this work. (Actually, 
it is empirically known that as long as the repulsion is sufficiently strong to avoid aggregation, the 
PCN model shows almost the same statistics.) 

To specify the PCN dynamics completely, we also need the network reconstruction rules. 
Whether the network is reconstructed or not is determined by the number of segments in chain 
end bonds. If the number of segments in a chain end bond becomes larger than a certain criterion, 
a new pair of slip-linked nodes are constructed. On the other hand, if the number of the seg- 
ments in an chain end bond becomes smaller than another criterion, the slip-link attached to that 
bond is destructed. The criteria are given as 3No/2 and iVo/2, with Nq being the average num- 
ber of segments in a bond. Although the network reconstruction rules affect the statistics rather 
strongly [10], how they affect the statistics is not so clear. (It is reported that other reconstruction 
rules [inilll] can improve several statistical properties.) In this work we do not consider about 
the reconstruction process and concentrate on the dynamic equations for {Ri k} and {Ni k} (eqs 



2.2 Detailed Balance Condition in the PCN Model 

Because the PCN model is designed as a dynamical model, its static properties are not clear from 
its dynamic equations (HJ-Q. Moreover, the existence of the thermodynamic potential (the free 
energy) is generally not guaranteed for such a phenomenological dynamical model. If eqs ([2])- ([4]) 
satisfy the detailed balance condition, the forces are expressed as variations of the free energy. (In 
the foUowings, we call such forms as the variational forms.) Then we can construct the free energy 
from eqs ©-([l]). In this section we attempt to interpret the PCN dynamic equations as Langevin 
equations _18^ which satisfy the detailed balance condition. 

According to the standard nonequilibrium statistical physics p"3llT9ll20] , dynamics of the system 
in the vicinity of equilibrium can be described well by the Langevin equations. If we assume that 
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there is no dynamic coupling between {Ri^k} and {Ni^k}, nor the memory effect, the detailed- 
balanced Langevin equations can be formally described as follows. 

-^f- - - 2. QR-^ + ^bTI^^ - i,,,^^., + e,fc {t) (9) 

3,1 • J,l 

Here i.-'^.j ; and i^^.^ ; are the mobility matrices, which may depend on stochastic variables such as 
{Ri^k} or {Ni^k}- (From the Onsager's reciprocal relation, these mobility matrices are symmetric.) 
J-eft({-Ri,fc}, {Ni^k}, {Zi}) is the effective free energy of the system, and the (generalized) forces are 
expressed as the derivatives of Tcs- and £,l^\t) are the Gaussian random noises which 

satisfy the the fluctuation-dissipation relations of the second kind. They can be expressed as 
follows, 
by eqs 
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where represents the transposed matrix of B. In eqs ([9]) and (fTO|) . the stochastic terms 
are interpreted following the Ito calculus [TB". The second terms in the right hand sides of eqs 
([9]) and pO)) are spurious (noise-induced) drift terms, which are required to satisfy the detailed 
balance condition (for the Ito form Langevin equations) . The dynamics described by eqs ([9]) is (fTO|) 
guaranteed to reach the equilibrium state, of which statistics is simply determined by the effective 
free energy J^eff- 

The PCN dynamic equation for {Ri k} (eq ^) seems to be similar to the Langevin equation 
([9]) . Actually, it is straightforward to show that eq ([2]) can be reduced to eq ^ by setting 

LlS3,i = ^^^^^^'^ (13) 
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and 
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Here J^cff,o({-Ri,fc}, {^i,fc}, {^i}) is the free energy for ideal linear springs, and J"cff,i({A'^i,fc}, {Zi}) 
is the remaining contribution to the free energy which is independent of {Ri^k}- 

On the other hand, the PCN dynamic equation for {Ni^k} (eq ([3]) together with eq (U))) is 
different from the standard Langevin equation PH)) . The segment number flux in the PCN dynamic 
equation, eq (U), can be rewritten as follows. 



Ji,k[t) — Mi,k 
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Substituting eq (jl7p into ([3]), we have the following expression for the dynamic equation. (For 
simplicity, we consider the case of 2 < /c < — 1 in eq ([3|). The generalization to fc = 1, is 
straightforward . ) 
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Comparing the stochastic terms in eqs (IT51) and (ITU)) , we find that the mobility matrix L^'^^j ^ 
be expressed by using Mi^k as 
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We also find that the deterministic terms in eq (jl9p can not be reduced to the variational form in 
eq ([TOl) . Namely, the PCN model does not satisfy the detailed balance condition and thus it does 
not have the thermal equilibrium state. The steady state in the PCN model (in absence of external 
flow field) corresponds to a nonequilibrium steady state, where the energy input and the energy 
dissipation is balanced pp] . 

Nevertheless, PCN simulations by the original dynamic equations can reproduce some equilib- 
rium properties of entangled polymers reasonably [1^. This implies that, although the PCN model 
is not detailed-balanced and does not have the equilibrium state, it works as a good approximation 
of a statistical mechanically sound, detailed-balanced dynamical model. Keeping this in mind, we 
consider to remedy the PCN model to satisfy the detailed balance condition in a heuristic way. 
From eqs ^T9\\ and ([10]), the variational form for the force can be realized if the following relation 
holds. 
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Clearly, eq (|2T|) does not hold if we use eq ([5]) as the definition of pi^k- Conversely, if we do not use 
eq ([5]) as the definition of pi^k, it becomes possible to satisfy eq ([21]) . Besides, there is no special 
reason to employ the arithmetic average form for pi_k- (Different models such as a constant value 
model were also proposed and utilized [1 .) Thus, here we define the segment density along the 
chain pi^k to recover the variational form. 
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The right hand side of eq (j22l) is the harmonic average of Ni_kJril\Ri,k+i — Ri,k\ and Ni^k/\Ri,k — 
Ri^k-i\, whereas the right hand side of eq ([5]) is the arithmetic average of them. By using the new 
definition ([22]) . it is straightforward to show that the condition ([21]) is satisfied and TcS^i = 0. 
Now the PCN dynamic equation for {Ni^k} (eq ([T9|) ) can be rewritten as follows. 
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Although eq ([^ is similar to eq ([TU)) . the spurious drift term (the second term in the right hand 
side of eq ([10]) ) is missing in eq ([23]) . Therefore, in addition to employ a new definition of pi^k (eq 
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([22])), we add the spurious drift term to the segment number flux equation (eq ([4])). We modify the 
definition for the segment flux (eq (|4])) as 



\Ri,k+l — Ri,k\ |-Ri,fc — -R'i,fe-1 I 



+ (25) 



where v\'^ is the spurious drift velocity defined as 
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The modified PCN dynamic equation for {Ni,k} (eq ^ together with eqs ([22]), ((251), and (|26| ) 
reduces to eq PH)) . Thus the modified PCN dynamic equations satisfy the detailed balance con- 
dition. It should be emphasized here that the noise term in eq (|IOp is the multiplicative noise 
because Lfl^j i depends on {Nij^} (through Mi ^ and Pi^k)- Although the spurious drift term may 
not be intuitive, it naturally arises as a property of the multiplicative noise |13[|21| , and is required 
to satisfy the detailed balance condition. 

Before we proceed to detailed analysis, here we shortly comment on the relation of our model 
to the GENERIC framework in the nonequilibrium thermodynamics [22ll23j . In our modified PCN 
dynamic equations, the thermodynamic forces are coupled to the symmetric mobility matrices. 
Such dynamic equations can be straightforwardly mapped onto the GENERIC framework. The 
mobility matrices and the effective free energy in our model correspond to the friction matrices 
and the entropy (with negative sign), respectively. They form the irreversible, dissipative bracket 
parts. The reversible, Poisson bracket parts do not exist in the modified PCN dynamic equations. 
Most of important properties in the nonequilibrium thermodynamics, such as the positivity of 
the entropy production, are automatically reproduced in the modified PCN dynamic equations. 
(Notice that, however, this does not mean that the modified PCN model is fully consistent with 
the GENERIC framework and the nonequilibrium thermodynamics. The repulsive interaction or 
the network reconstruction process may be inconsistent with them.) 



2.3 Effective Free Energy 

The modified PCN dynamic equations introduced in the previous subsection reduces to the detailed- 
balanced Langevin equation. The equilibrium state of the modified PCN model is determined by 
the effective free energy J^^s- The effective free energy of the modified PCN model can be simply 
expressed as 

TMiR^.k}, {N^.k}, {Z.}) = ^eff,o({i2..}, {N.,k}, {Z^}) = ^ E E ^^J'"''^' (27) 

i k=l ' 

Notice that, the effective free energy ([27)) is different from the free energy of ideal Gaussian chains. 
The free energy of Gaussian chains is expressed as 

M{R.kh {N.^kh {Z.}) - ^ E E ^""^ 'n^^' + ^ E E '^-^>^ (28) 

i k—1 i k—1 

Comparing eqs (f27| and ([28]) . we find that the logarithm term is missing in the PCN effective free 
energy. 

We can rewrite the PCN effective free energy (|27p as follows, by using the free energy of 
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Gaussian chains ((28)) . 

Zi 

T,S{{R^A^ {Zi\) = :F^{{R^,k}, {N,,k}, {Z,}) + ^ sMpA,n\.{N k) (29) 

i k=l 

C^siip-iink(?i) = - In n (30) 

Here Usii-p-\i-ak{n) represents the effective interaction between neighboring shp-hnks. (Strictly speak- 
ing, L'siip-iink also represents the effective interaction potential between a slip-link and a chain end. 
In the foUowings, we assume that the interaction potentials between the slip-links and between a 
slip-link and a chain end are always the same.) We note that ?7siip-iink("-) depends only on the 
number of segments between slip-links. This means that the effective interaction between slip-links 
is determined by the chemical distance, not by the spatial distance. Eq ([50]) is monotonically de- 
creasing function of n and thus the interaction between slip-links is purely repulsive. In addition, 
the effective potential diverges at the limit of n — 0, and thus slip-links strongly repel each other 
at the short range (at the short chemical distance). 

The second term in the right hand side of eq ([29)1 can be interpreted as the total repulsive 
interaction energy of slip-links. In this sense, we can say that the PCN model is based on the 
Gaussian chain and repulsive slip-links. It is worth mentioning that such an effective repulsion 
between slip-links is already suggested by the CReTA primitive path extraction for Monte Carlo 
simulation data [5]. 



3 Equilibrium Probability Distributions by Single Chain 
Approximation 

3.1 Single Chain Approximation 

In the previous section, the effective free energy for the modified PCN model (l?7|) is shown to be 
different from the free energy for Gaussian chains (j28l) . due to the repulsive interaction between 
slip-links. This means that several equilibrium statistical properties of the PCN model are qualita- 
tively different from ones of the Gaussian chains with slip-links without any direction interactions. 
However, unfortunately we cannot calculate the equilibrium statistics of PCN analytically from 
the effective free energy. (It is practically impossible to calculate the statistical weight for strongly 
slip-linked multi chain systems.) In this section, we calculate the equilibrium statistical properties 
by using the mean field type single chain approximation. (We notice that, it is not clear whether 
the single chain approximation employed in this work is really reasonable for the PCN model 
or not. We employ the simplest approximation to make the expressions simple and analytically 
tractable.) 

Schieber [241 proposed a single chain model with slip-links in which slip-links behave as a sort 
of grand canonical gas particles on a polymer chain. In his model, the number of slip-linked 
subchains is controlled by the effective chemical potential. There is no direct interaction between 
slip-links and the polymer chain obeys the Gaussian statistics. Namely, slip-links are placed on 
the target polymer chain and they do not directly interact each other via a potential (this model 
corresponds to a sort of ideal model). In this work, we call his model as "the single chain non- 
interacting slip-link model" or simply as "the non-interacting slip-link model" . (The expression 
"non-interacting" may sound somehow unnatural, because the slip-links constrain the polymer 
chain and this constraint effect can be interpreted as a sort of interaction. Nonetheless, in this 
work we employ this expression because in the followings we consider another model in which 
slip-links directly "interact" each other.) The statistics of the single chain non-interacting slip-link 
model is briefly summarized in Appendix [X] 

Here we follow Schieber's idea and consider the single chain version of the modified PCN 
model. We may call this as the "the single chain repulsive slip-link model" or "the repulsive slip- 
link model" in the followings. This is because in our model, neighboring slip-links interact each 
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other via the repulsive potential (eq ([30]) V In the single chain repulsive slip- link model, the state of 
the system is expressed by the node positions {Rk} {k represents the A:-th segment on the chain), 
the numbers of segments {Nk}-, and the number of subchains Z . (Under the mean field type single 
chain approximation, the connectivity map is not required.) 

The effective free energy for a single Gaussian chain with repulsive slip-links can be expressed 
as the single chain version of eq ((27| . 



^eff(TO,™,^) = ^E^ 



R 



■k-lj 



fe=l 



(31) 



Introducing the effective chemical potential for subchains, e, we define the following effective grand 
potential (grand canonical thermodynamic potential) for the repulsive slip-link model. 



JesiiRk}, {Nk}, Z) = F,s{{Rk}, {Nk}, Z) - eZ 



(32) 



The effective chemical potential e is determined from the following condition for the equilibrium 
average number of subchains 

N 



{Z) 



cq 



(33) 



where 



/cq 



represents the equilibrium statistical average calculated by the equilibrium prob- 



ability distribution. Here N is the number of segments in a chain and TVq corresponds to the 
equilibrium average number of segments in a subchain. 

All the equilibrium statistical properties can be calculated by the effective grand potential (|32p . 
The equilibrium statistical probability distribution is given as the standard Boltzmann distribution. 
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(34) 



where A is the thermal de Broglie wave length and S is the grand partition function. (The subscript 
"eq" in eq (j34p refers to the equilibrium state.) The thermal de Broglie wave length was introduced 
to make the grand partition function dimensionless. The delta function in eq (j34p comes from the 
constraint that the total number of segments in a chain is constant. The grand partition function 
S is defined as 



^^T^jmrk^ f d{Rk}d{Nk}s(N~Y.^^)^^p 

z=\ •' ^ fe=l ^ 



J,s{{Rk},{Nk},Z) 



Here we introduced the shorthand notation for the integrals over {Rk} and {Nk}- 



z „ Z „7V 
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By calculating the integrals over {Rk}, eq ((35)) becomes as follows. 
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where V is the volume of the system. The integrals over {Nk} in eq (|37l) can be analytically 
calculated as 



ATSZ/ 



\^{d{Nk}dN-Y^N^Y^N\ 

^ fc=i ^ fe=i 
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(38) 



where T{x) is the gamma function [25]. The detailed calculation is described in Appendix lB.il By 
substituting eq (|38l) into eq (|37l) . finally we have the following expression. 
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Here we defined the dimensionless effective fugacity ^ as 

^ is determined to satisfy the condition for {Z)cq (eq ([551) 1. Ea.pix) in the last hnc of cq is 
the (generahzed) Mittag-Leffier function [26] defined as 

n— 

Unfortunately, the Mittag-Lefflcr function is not easy to handle. Therefore, in the following sub- 
sections we derive simple and approximate expressions rather than exact expressions. 

3.2 Slip-Linked Subchain Number Distribution 

First we calculate the slip-linked subchain number distribution. The distribution function of the 
subchain number Z can be expressed as follows. 

Peq(^) = I d{Rk}d{Nk} Pcq({i?fc}, {Nk}, Z) (42) 

By integrating eq (|42|) over {Rk} and {-/Vfc}, we have the following expression. 



Eq (|43|) can not be expressed by elementary functions. This is in contrast to the non-interacting 
slip-link model, of which subchain number distribution is simply given as a Poisson distribution 
(eq (US]) in Appendix El . 

We can utilize the saddle point approximation to obtain a simple approximate form for eq (j43p . 
(See Appendix [Un i 

y2—{Z - Zq) 



l5Zo + l 



4Z^ 



(44) 



Figure [T] shows the subchain number distribution functions for various values of Zq by eq (j43p 
(exact) and by eq (^1)) (approximation). For comparison, the distribution function of the non- 
interaction slip-link model (eq ((75)) in Appendix |X| is also shown in Figure [TJ We find that the 
approximate form (j44p works very well even for relatively small Zq such as Zq = 5. 
The relation between the average and variance of Z is easily obtained from eq (|44)) . 



Oy2 o 

{{Z-Z,fU^^^-^^-Zo (forZo»l) (45) 

Figure [2] shows the relation between {Z)eq = Zq and {{Z — Zo)^)eq, calculated by the exact and 
approximate forms. Again we find that the approximate form ()44p works well even for relatively 
small Zq. For sufficiently large Zp, we find that the variance of the repulsive slip-link model is 
smaller than one for non-interacting slip-link model by the factor 2/5 (see eq (|74l) in Appendix [5| . 
The distribution of Z of the repulsive slip-link model is sharper than one of the non-interacting 
slip- link model, as shown in Figure [T] 

3.3 Segment Number Distribution 

Second we calculate the distribution function of the segment number in a subchain. We express 
the number of segments as n. The segment number distribution function Poq(n) is expressed as 
follows. 



Poq{n) ^Y. d{Rk}d{Nk} 
z=i •' 



z 

-Y^6{n-Ni) 

1=1 



P,q{{Rk}ANk},Z) (46) 
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By calculating the integrals over {Rk} in eq (|46l) . we have 



V 



E 



f [I E -^(^ - A^ol -5 (a^ - E ^0 n ^- 

Lii J\ I — 1 / I — 1 



d{Nk} 

' 1 = 1 -' ^ k=l ' k=l 

Integrals over {Nk} in eq (|T7)) can be calculated analytically. For Z > 2 we have 

z n / z z 



3/2 



(47) 
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,ix5(z-i)/2-i (30F/4) 



AT 



(48) 



r(5(Z-l)/2) 



(see Appendix IB.ip and for Z = 1 we have 



1 



7V1/2 

Thus eq (|47p can be rewritten as follows 
4 



N 

dNi S{n - Ni)d{N ~ Ni)Nf^^ ^ ^ {'^ 



^(^-)-i:(^)'"(-^) 



5(Z-l)/2-l 



Z-1 



Z=2 



r(5(Z-l)/2) 



(49) 



(50) 



As before, we can obtain a simple approximate expression by using the saddle point approxi- 
mation. (See Appendix IC. 21 ') For sufficiently large Zq, eq (|^5)) is approximated as 



25 /l0n3/2 



6 



exp 



5n 
'2]Vo 



(51) 



Eq (jSip is qualitatively different from one of the single chain non-interacting slip-link model (eq 
(j78p in Appendix E| . Figure [3] shows the segment number distribution functions for repulsive and 
non-interacting slip-link models, for sufficiently large Zq. The distribution functions are normalized 
to satisfy the following normalization condition. 



din/No) Peq{n/No) = 1 



(52) 



We can observe differences between two models clearly. The difference is especially large for small 
n. Eq (I5ip is not a monotonically decreasing function and it has a maximum at u/Nq = 3/5. Also, 
eq ([5T|) approaches to at the limit of n — ^ 0, while the distribution function of the non-interacting 
slip-link model approaches to a non-zero constant. Intuitively, this can be understood as the effect 
of the strong short range repulsion between two neighboring slip-links. 



3.4 Bond Vector and Bond Length Distributions 

Finally we calculate the distribution functions for the bond vector and bond length. Here we define 
the bond vector as the vector which connects the two successive slip-links (or pair of a slip-link and 
a chain end) on the same chain. If we express the bond vector as Q, the bond vector distribution 
function Poq(Q) is expressed as follows. 



PoM = E / d{Rk}d{Nk} -Y,6{Q-Ri+ Ri I 

z=i'' L 1=1 



P,^{{Rk}ANk},Z) 



Eq ([55)) can be modified as follows. 
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The integrals over {Nk} can be calculated in the same way as the case of the segment number 
distribution function. Then we have 

For sufficiently large Zq, the bond vector distribution function can be approximately expressed as 
follows, by using eq ([5T|) . 



(56) 



6 V TT V 2TrNob^ J Nq Jq \ 2nb^ 2No 



2tt^ N^b-i \\ Nob^ 



where Ki{x) is the first order modified Bessel function of the second kind [25 . (The detailed 
calculation is described in Appendix IB.21 ) From eq (|56p . the bond length distribution function 
(the bond length Q is defined as Q = \Q\) becomes 



PM ^ ^.Q^PM - ^l^i^i {^q) (57) 

From eq ()57p we find that the form of the bond length distribution function is also different from 
one of the single chain non-interacting slip- link model (eq (|8ip in Appendix]^. Figure 2] shows the 
bond length distribution functions for repulsive and non-interacting slip-link models, for sufficiently 
large Zq. The distribution functions are normalized to satisfy the following normalization condition 



diQ/^N^) Pcq(g/ = 1 (58) 



For small Q, the asymptotic form of eq ([57|) becomes 

10Vl5 Q 



-Pcq(Q) « — fn- « (for Q « VnJ^) (59) 

TT N^' b^ 

This dependence of the asymptotic form on Q is different from one for the non-interacting slip-link 
model (PcqiQ) oc Q). The difference between two distributions is especially large for small Q. This 
trend is similar to the difference between the segment number distribution functions. The position 
of the maximum of eq dST]) is Q/y/N^ sa 0.616 whereas it is Q/y/NJ^ = 1/VG ~ 0.408 for the 
non-interacting slip-link model. As before, these properties can be understood as the effect of the 
strong short range repulsion between slip-links. 



3.5 Comparison with PCN Simulations 

The distribution functions obtained in preceding subsections (Sections 13. 2113. 4p are based on the 
single subchain approximation. However, the validity of the single chain approximation is not 
guaranteed for the PCN model. In this subsection, we perform PCN simulations with original and 
modified dynamic equations and directly calculate the distribution functions. 

The simulations are performed for linear monodisperse polymers with Zq = 5, 10, 20 and 40. 
The osmotic terms to prevent aggregation of slip-links [31[TU] are added to eqs ©, (H)), and PT]! . 
The network reconstruction is used to realize the steady state (that is the equilibrium state for the 
modified model). The employed rule was NRl in Ref [TOl fwhich is briefly explained in Section [OJ . 
The unit cell dimension for periodic boundary condition is and the total bond number 

in the cell is ca. 10 x 12'^. The distribution functions are obtained for 10 independent simulations 
and the averages are reported here for noise reduction purposes. 

Figure [5] shows the comparison of the subchain number distribution functions Peq{Z) calculated 
by the single chain models and the PCN simulations. As shown in Figure [SJ the data of the 
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modified PCN model agree well with the repulsive slip-link model. The original PCN model gives 
a subchain number distribution function which is broader than the repulsive slip-link model but 
still sharper than the non-interacting slip-link model. Considering the roughness of the single chain 
approximation, the agreement between the modified PCN model and the repulsive slip-link model 
is rather surprising. The CReTA primitive path extraction by Tzoumanekas and Theodorou [3] 
also shows subchain number distributions which are sharper than the Poissonian, but still broader 
than eq (|44)) . Thus these models will be in between the non- interacting slip-link model and the 
repulsive slip-link model (if we assume that the repulsive interaction between slip-links is the main 
reason which cause non-Poissonian distribution). 

Figure ini shows the comparisons of the segment number distribution functions Peq{n) for Zq = 
40. (The sharp peaks of Pcq{n) by the PCN simulations around n/No = 0.1 are the artifact due 
to the cutoff of the segment number JlO^. In the foUowings we neglect these artificial peaks.) 
Unlike the case of the subchain number distribution, we can observe the deviation of the modified 
PCN model data from the repulsive slip-link model. As shown in Figure [6l the segment number 
distribution function of the modified PCN model is somehow shifted to large n region. However, 
the asymptotic behavior of the modified PCN model at small n is consistent with the repulsive 
slip-link model {Peq{n) oc n^/^). The original PCN model gives the distribution between the non- 
interacting and the repulsive slip-link models. This seems to be qualitatively similar to the case 
of the subchain number distribution. Here it is worth mentioning that the CReTA primitive path 
extraction data [9] gives the qualitatively similar asymptotic behavior, P{n) — >■ at n — 0, to 
the modified PCN model and the repulsive slip-link model. We will discuss the segment number 
distribution functions later in detail (Section 14. 2[) . 

Figure [7] shows the bond length distribution functions Pcq(Q) for Zq — 40. We can find the 
bond length distribution function of the modified PCN model deviates from the repulsive slip-link 
model. This result is physically natural, because the bond length distribution function and the 
segment number distribution function are expected to be related each other. In the single chain 
model, actually their relation is simply expressed by eq (j55l) . The bond length distribution function 
of the modified PCN model is slightly shifted to large Q region. The asymptotic behavior of the 
modified PCN model at small Q is also consistent with the repulsive slip- link model (PcqiQ) oc Q^). 
The original PCN model data is again between the non-interacting and repulsive slip-link models. 

We can observe that both the original and modified PCN simulation data have relatively large 
probability for large Q. This implies that subchains in the system are stretched compared with the 
ideal Gaussian subchains. This would be because chains are required to form a connected, tetra- 
functional network structure in the PCN model. (In the PCN model, a slip-link spatially binds two 
chains, and the bound chains will be somehow stretched.) This can be interpreted as a result of 
the interchain force-balance condition imposed in the PCN model (this condition is automatically 
imposed by eqs ([T]) and ©). The possible chain conformations are limited compared with the 
single chain model, where no such constraint exists. As a result, the bond length distribution 
function becomes broader and deviates from the distribution function in the single chain model. 
This is qualitatively consistent with the simulation result. Also, the segment number distribution 
would shift to large segment number direction, which is again qualitatively consistent with the 
simulation result. Therefore, we expect the constraint for chains to be the origin of the deviation 
of PcqiQ)- The origin of the deviation of Pcq{n) will be also the same. 

In the single chain model, effects of the surrounding chains are expressed only by the effective 
chemical potential e. Such a rough approximation would not be fully applicable to describe the 
statistics of slip-linked multi chains. The effective chemical potential controls only the number of 
slip-linked subchains, and the statistics of bond vectors or segment numbers are determined only by 
the effective free energy. If interactions between the target chain and the surrounding chains are not 
simple and cannot be described by a single parameter e, the single chain model cannot reproduce 
the PCN simulation results. To analyze the statistical properties of the PCN model theoretically, 
therefore, we will need to improve the single chain approximation to take account the effects of 
surrounding chains. For example, some sort of spatial correlation effect may be took into account 
to the model, by using the connectivity information (such as the connectivity map C in the PCN 
model). Several radial distribution functions of nodes are obtained by the CReTA primitive path 
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extraction [9]. Kindt and Briels [21] modeled that the average number of entanglements between 
two chains as a function of the distance between centers of mass of chains for their highly coarse- 
grained single particle model of entangled polymers. If we use of the connectivity information and 
reproduce these spatial correlations, several statistical properties may be improved. Of course, 
there are other possible factors which affect the statistical properties of the PCN model. (As we 
mentioned, the network reconstruction rule [10] affects several statistical properties of the PCN 
model. Recent simulation results show that the spatial correlation between nodes is somewhat 
affected by the form of the interaction potential between nodes, although rheological properties 
are almost unchanged [27]. Such a local structural change may affect the distribution functions 
to some extent. If we directly control the spatial correlations of slip- links by controlling the 
reconstruction rules and the connectivity information, the distribution functions will be affected.) 



4 Discussion 



4.1 Detailed Balance Condition and Effective Free Energy 



As shown in Section 12.21 the PCN dynamic equations can be reduced to the detailed-balanced 
Langevin equations with two modifications. It is rather surprising that simple modifications 
can reproduce rather complicated, detailed-balanced Langevin equations. The modified (detailed- 
balanced) version of the PCN model can reproduce several probability distribution functions which 
are qualitatively similar to the ones predicted by the single chain repulsive slip-link model. 

Here we consider the original PCN dynamic equation for {Ni,k} (eqs (O and (dj)) again. The 
fact that we only need relatively minor modifications for the original PCN dynamics equations 
implies that the original PCN dynamics itself is already a reasonable approximation for a statistical 
mechanically sound model. Actually, we can show that the original PCN model approximately 
satisfies the detailed balance condition under a certain condition. We consider the most probable 
value of Ni^k under a given {Ri^k} and Zi. The most probable state {A"i,fc} minimizes the effective 
free energy (|?7|) . Ni^k is simply given as 



Ni^k — Pi\RiM — Ri,k-l\ 



(60) 



where pi is the most probable value of the monomer density along the chain on the i-th chain {pi^k 
is independent of k and thus we simply describe it as pi). If the deviation of Ni^k from the most 
probable value is small, we can rewrite Ni^}^ as TVi.fc — Ni,k + SNi,k and expand the free energy or 
the dynamic equation into power series of {SNi^k}- Expanding eqs ([^ . ([^ and ([?f)) into power 
series and retaining only the leading order terms, we have 
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(61) 



(62) 



(63) 



The original PCN model, eqs ([5]) and (|3|) also reduce to eqs ([5^ and ([M]). respectively. Therefore, 
the original PCN model satisfies the detailed balance condition if {Ni^k} is in the vicinity of the 
most probable state {Ni^k}- (Intuitively, this is because the segment density pi^k is not sensitive to 
the node index k in the vicinity of the most probable state. The harmonic and arithmetic averages 
are almost the same and the spurious drift velocity is negligibly small, under this condition.) 
Unfortunately, this condition is not always satisfied even in equilibrium, because SNi ^ is generally 
not negligibly small due to the thermal fluctuation. Nevertheless, the probability that {Nij^} is 
close to {Ni^k} is expected not to be small. Thus we consider that the analysis performed here 
qualitatively explains why the original PCN model works as a good approximation. 
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The new definition for the segment density pi^k (eq (|22|) ) proposed in this work is expressed 
as the harmonic average of the local densities in two neighboring bonds. Although this form is 
not intuitive, from the view point of the discretization scheme, it is not unreasonable. The PCN 
dynamic equation for {Nt^k} (eqs ([3|) and (j4|)) is interpreted as the stochastic diffusion equation 
on one dimensional discrete lattice. At the continuum limit, the PCN dynamic equation will 
reduce to the stochastic partial differential equation |28| which contains the multiplicative noise. 
To discretize the multiplicative noises appear in the stochastic partial differential equations, non- 
arithmetic averages are sometimes preferred than the simple arithmetic average. For example, the 
multiplicative noise in the dynamic density functional equation can be accurately discretized with 
the geometric average form [33] . (The geometric or harmonic average is zero if one of the averaged 
numbers is zero, whereas the arithmetic average is nonzero. Such a property is required to generate 
the multiplicative noise accurately and stably.) Thus we expect that the harmonic average form is 
more natural than the arithmetic average form. 

We should note that the statistical properties of the PCN network reconstruction rules flOjlTT] 
are not considered in this work. The network reconstruction rules proposed so far do not satisfy 
the detailed balance condition and therefore the PCN model is not fully detailed-balanced, even if 
we employ the modified dynamic equations. (Especially for small Zq, the network reconstruction 
rule affects the statistical properties significantly.) To understand the statistical properties of 
the full PCN dynamics, the statistical mechanical analysis of the network reconstruction rules 
is also required. It seems to be difficult to satisfy the detailed balance condition in the network 
reconstruction process of multi chain slip-link models due to its intrinsic complicatedness, and thus 
it is left for a future work. (The exact equilibrium statistics of a slip-linked multi chain system is 
not simple nor intuitive |16) . unlike a single chain system.) Although we do not consider further 
modifications here, the current work will be useful to find strategies for other modifications which 
improve the PCN model. 

4.2 Segment Number Distribution 

As we showed and discussed in previous sections, the effective repulsive interaction between slip- 
links (eq (PH)) ) affects statistical properties of the network qualitatively. We expect that the statis- 
tics can be further changed by changing the interaction potential (or interaction strength) between 
slip-links. In fact, if we consider the limit of the strong repulsion, distribution functions become 
much different from non-interacting slip-link model (see Appendix ID|) . Roughly speaking, the 
variance of the subchain number distribution function can be used to estimate the strength the 
effective repulsive interaction between slip-links. It is interesting to compare our model (or our 
picture) with other models or primitive path extraction results, and discuss whether our model is 
similar to others or not. 

Among the examined distribution functions, the segment number distribution function seems to 
be the most sensitive to the interaction between slip-links. The segment number distribution of the 
single chain repulsive slip-link model is qualitatively much different from one of the single chain 
non-interacting slip-link model. The segment number distribution functions which have similar 
forms to eq ((5T|) are already proposed by Tzoumanekas and Theodorou [9], and Greco [30] . 

Here we compare our model with these distribution functions. Tzoumanekas and Theodorou 
proposed the following empirical form for the segment number distribution, to fit their CReTA 
primitive path extraction data. 

^eq(n) = l^-^(e-^"/^° - e-'^"/^") (64) 

JVO c — 

Here 1 < 6 < 2 is a fitting parameter and c — b/{b—l). They reported that the segment number 
distribution functions obtained by the CReTA can be fitted well to eq (|64| with h — 1.30. Eq (|64]) 
can be also utilized to fit the distribution functions by the original PCN model with the fitting 
parameter b — 1.96 10 . This value is close to 2, and for such a case eq ([M| can be approximated 
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well by the following simple form. 



An 



-2n/Na 



(65) 



Eq ([55]) is quite similar to eq (|5ip . Tzoumanekas and Theodorou proposed that the physical 
origin of such a non-monotonic form of eq (j64p is the effective repulsion (or the blocking effect) 
between topological constraints, which depends the chemical distance. This is qualitatively the 
same as our repulsive slip-link model. Of course, the CReTA primitive path extraction results and 
their molecular model are different from the PCN model and our slip-link model. For example, 
the CReTA primitive path data are taken for an atomistic polymer model while the PCN model 
employs a rather coarse-grained Gaussian chain model. Besides, it is not clear whether the slip- 
links in the PCN model and the topological constraints extracted by the CReTA have the same 
physical properties. Thus the distribution function by the PCN model does not necessarily to be 
coincide to one by the CReTA. Actually, at least the values of the fitting parameter b for two 
distributions are different. Nonetheless we consider that the similarity between the Tzoumanekas- 
Theodorou model and our repulsive slip-link model implies that there are some common properties 
of slip-links (or topological constraints). (We expect that the detailed form of the effective repulsive 
potential between topological constraints in the Tzoumanekas- Theodorou model is different from 
our repulsive slip-link model. Judging from the distribution functions, the repulsive interaction 
in their model seems to be weaker than our model. But the detail repulsion mechanisms in the 
Tzoumanekas-Theodorou model or the CReTA results are not so clear. Further investigations will 
be future works.) 

On the other hand, in the Greco model, the segment number distribution function is given as 



Peq(n) oc si?n)\j—e 
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(66) 



where m > 1 is the renormalized chemical potential (which is used as a fitting parameter), s{x) = 
(x^ — i)/x, and erfc(a;) is the complementary error function 25 . Eq (|66p has a similar form to eq 
([64)) |10P30]. In the Greco model, the slip-links are non-interacting. The non-monotonic nature of 
eq comes from the fiuctuation effect due to the smallness of the system. This is qualitatively 
different from our model. It is interesting that different mechanisms (the fluctuation effect and the 
repulsive interaction between slip-links) give the similar result. Further progress of theories would 
be required to understand the origin of non-monotonic segment number distribution functions. 

Although there will be other factors (the osmotic force between nodes, the network recon- 
struction, or constraint to form the tetra- functional network structure), the interaction between 
slip-links is one important factor which determines the network statistics. The concept of the 
effective repulsive interaction between slip-links will be also useful to consider statistical properties 
of other models or methods such as the primitive path extraction methods f7}^ or other slip-link 
based models without explicit free energy model |31| . For example, the differences among various 
primitive path extraction methods may be related to the differences of the effective slip-link inter- 
action potential models. Once we obtain the explicit form of the effective interaction potential, it 
will be possible to tune it so that the model reproduces the required equilibrium statistics. We can 
employ the effective interaction potential for slip-links determined from the network statistics by 
the primitive path extraction methods to make the slip-link models compatible with the primitive 
path extraction data. (We can easily replace the effective free energy model by other models to 
tune the statistics of the modified PCN model, because the forces in the modified PCN model are 
expressed as variational forms.) 



5 Conclusions 

We showed that the PCN dynamic equations do not satisfy the detailed-balance condition and 
therefore the PCN model does not have the thermal equilibrium state. However, by introducing 
heuristic modifications, the PCN dynamics can be recovered to satisfy the detailed balance condi- 
tion. We proposed two modifications. One is to change the definition of the segment density along 
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the chain to the harmonic average form (eq (|22]) ). and another is to add the spurious drift term to 
the monomer flux equation (eq ([25|) ). From the modified PCN model, we obtained the effective free 
energy from which the equihbrium statistical properties of the PCN model are determined. The 
effective PCN free energy has a different form from the free energy of ideal Gaussian chains. This 
can be understood that the effective PCN free energy contains the contribution from the repulsive 
interaction between slip-links. 

To analyze the equilibrium statistical properties of the PCN model, we constructed a single 
chain model with repulsive slip-links. The equilibrium distribution functions derived from the 
single chain repulsive slip-link model are qualitatively different from ones for the single chain 
non-interacting slip-link model, due to the repulsive interaction between slip-links. It was shown 
that the slip-linked subchain number distribution function for the repulsive slip-link model is 
much sharper than one for non-interacting slip-links. The interaction between slip-links strongly 
affects the segment number distribution function, especially if the segment number is small. The 
equilibrium properties obtained from the single chain slip-link model are qualitatively similar to 
the PCN simulation results, although the agreement is not perfect. We consider the repulsive 
interaction between slip-links is an important factor which determines the statistical properties of 
the PCN model. The repulsive slip-link picture would be useful to improve the model statistics or 
design a new model. 
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Appendix 

A Single Chain Non-Interacting Slip-Link Model 

In this appendix, we show a brief derivation of the equilibrium distribution functions for the 
single chain non-interacting slip-link model [24] . (The derivation mainly follows Ref [24| but the 
expressions are slightly modified to allow the direct comparison with the repulsive slip-link model 
in the main text.) The free energy of a Gaussian chain with non-interacting slip-links is expressed 
as the single chain version of eq 
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As we mentioned in the main text, the Gaussian free energy contains a term which is proportional 
to IniVi while the PCN effective free energy ([ST)) does not. 

The grand potential is defined as follows, by introducing the effective chemical potential for a 
subchain, eg. 

M{Rk}, {Nk}. Z) EE Fa{{Rk}, {Nk}, Z) - e^Z (68) 



The equilibrium statistical probability can be written as 
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with ^0 being the grand partition function. The grand partition function reduces to a simple form. 
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Here we defined the dimensionless efi^ective fugacity as 



The slip-Unked subchain number distribution function becomes the Poisson distribution. 

Po^{Z) = J d{Rk}d{m} P,^i{R,},{N,},Z) ^ j^^^^t'e-^" 
Eq ((72)) can be rewritten by using the relation, (Z)cq = = Co + 1- 

The variance of Z is related to the average of Z as 
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The segment number distribution function is expressed as 
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The integrals over {Nk] can be calculated to be 
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Thus eq (|75p can be reduced to the following form. 
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where we used Zq — N/Nq (notice that, this Nq is different from in Ref [24|) . For larg 
(|77p can be approximated by the exponential function. 
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The bond vector distribution function can be calculated in a similar way. 

z 

P,^{{Rk},{Nk},Z) 



Pcq(Q) - E / d{Rk}d{Nk} -Y,^{Q-Ri + Ri I 

z=i-^ L 1=1 



27r62 



3/2 „Ar 



1=1 



17 



Substituting eq ([77)) into eq (|79)) . we have 
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where erfcx is the complementary error function f35]. (See Appendix IB. 21 for details.) Although 
eq (|80p is exact, it is quite complicated and not intuitive. It can be reduced to the following simple 
approximate form for sufficiently large Zo- 
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Nob^ 



IQI 



The bond length distribution function can be expressed as 

^cq(Q) 



Nob^ 



exp 



Nob^ 



Q 



(81) 



(82) 



Eq ([82]) has a maximum at Q/ \/NoF^ = 1/ \/6. 

As shown in the main text, the distribution functions for the single chain non-interacting slip- 
link model are qualitatively different from ones for the single chain repulsive slip-link model. Such 
differences arise only from the difference of the interaction between neighboring slip-links. 



B Calculations of Integrals 

In this appendix, we show the detailed calculations for several integrals appear during the derivation 
of distribution functions in the main text. Although the calculations are complicated, we can obtain 
the analytic expressions by using special functions such as the gamma function or the modified 
Bessel function. 



B.l Integrals over {Nk} 

First we show the detailed calculation of the integral which appears in eq ([38 

z . z 



^1 - ]^irT / d{N,}s(N-J2N,)Y[N!/' (83) 

^ fc=l ^ k=l 



To calculate the integral over {-/Vfc}, we introduce the following variable transform. 

fo (fc = 0) 

Sk = l{Nk- Nk-i)/N (fc = l,2,...,Z-l) (84) 
[l (k^Z) 

The integral over {Nk} with the delta function can be transformed into integral over {sk} without 
a delta function. 

dN,s(N-Y.^'^) = ]l (85) 
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Eq ([83]) can be modified as follows, by using the transform ((85 



fc=i 



z-i 



i=0 



We consider the integral over si in eq (|86l) . It can be modified as 

dsi (S2 - ' = 5(5/2, 5/2)4 

Jo 



(86) 



(87) 



where B{x, y) is the beta function Similarly, the integral over S2 in eq can be calculated 
to be 



/ ' ds2 (S3 - S2f^^ [5(5/2, 5/2)4] ^ B(5/2, 5/2)5(5/2, 10/2); 
Jo 

Iterating the same procedure, finally we have the following expression for Ii. 



13/2 



Z-1 



Z-1 



h = n^(5/2,5fc/2)= n 



r(5/2)r(5fc/2) _ r^(5/2) 
r(5(fc + i)/2) ^ r(5z/2) 



fc=l fc=l 
where wc have utilized the relation between the beta and gamma functions |25) . 

T{x)r{y) 



T{x + y) 



Bix,y) 



By substituting r(5/2) Sy/ir/A into eq we have eq (p8| . 
Next we show the detailed calculation of the integral in eq (jUl 



(89) 



(90) 



1 



Z Z-1 



r z 

•' L 1=1 



5{n-Y^nAX{n, 



k=l ' k=l 



1=1 k=l •'^ m=0 



(91) 



^7n + l -^m) 



3/2 



where we defined n — n/N. The integral over {sk} can be performed in the similar way to the 
case of Ii . 



z z-1 



-^2 = ^ XI n / cisfc n ~ *™)^^^ 

X dsi 5{fi - (sz+i - si)){si+, - sif'^sf-^^'^ n ^(W2, 5/2) (92) 

=^E n r^^'^^fc n (s™+i-5m)'/'(^m-r.)(5'-2)/2-Q 5(5^/2,5/2) 



1=1 k=l+l ' 



m=l+l 



m=l 



To simplify the expression, we introduce another variable transform, Uk ~ Sk — fi. Then we have 



^;3/2 



Z Z-1 

E n 

1=1 k=i+i 



z-1 l-l 

duk n - Umf'^4t;^'^'^ n B{hm/2, 5/2) 

7n— /+1 m— 1 

Z-2 



^x:(i--)''''^"^"''''n^(W2,5/2) 



(93) 



/=1 m=l 

Finally we have the following expression for 12- 

-^i3/2.-, _ -.n5(Z-1)/2-1_ 



^Z-1 



l2=n"\l-n) 



(5/2) 



r(5(Z-l)/2) 

Eq dM]) together with n = n/iV and r(5/3) = gives eq (|48| . 



(94) 
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B.2 Integrals over n 

Here we show the detailed calculation of eq ([56]). The integral over n in eq ([56]) can be modified 



as follows, by using the variable transform t = y^bb'^ /SNqQ^u. 



dn exp 



Jo 



3Q2 5 n 



3g2 



5Nob^ Jo 



dt exp 



1 /I5g2 /I 



iVo52 \t 



+ t 



(95) 



The integral in eq ((95)) can be reduced to the first order modified Bessel function of the second kind. 
The integral expression of the first order modified Bessel function of the second kind becomes [25] 



Kiix) 



due-^'=°"^"coshM = 



f 



X cosh u 



1 



POO 




1 ds exp 




/o 





(96) 



where we have used the variable transform s — e By substituting x ~ 15Q^ / NqB'^ into eq 
(|96p . we have the following expression for I3. 



h = -1 



'f5Q2 



Ki 



'f5 



Noby 



(97) 



This gives eq (l56t . 

A similar integral also appears in eq (IMl) . The integral in eq ([501 can be calculated as follows. 



N~ N, 



N 



NoN 7o 
2{Zo - 1) 

„3/2 



dn ( — ^ 



3/2 



dt exp 



exp 

■)2 



3Q2 N-No 
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NoN 
1 ' 



^^t'-iZo-l)^ 



(98) 



where we used the variable transform t = y/N/n. The integral over t in eq (|98]) can be calculated 
by using the following formula |25] . 



erfc(Va? + %/y) + e-2v^erfc(VJ - ^) 



(99) 



where erfc x is the complementary error function [25j . Finally eq ([98p can be reduced to 

erfc 



I4 = 



V^(Zo-l) Nb^ 



/3/2 



■ exp 
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(100) 



This gives eq (|80]) . 



C Calculations of Saddle Point Approximations 
C.l Saddle Point Approximation for PcqiZ) 

A simple approximate expression is demanding to analyze the statistical properties of the repulsive 
slip-link model. Here we attempt to calculate the approximate expression for eq ([43p , for sufficiently 
large Zo- For convenience, we approximate Pcq{Z) by a continuum distribution. The summation 
over Z is replaced by the integral as follows. 

^ « / dZ (101) 
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We can utilize the Stirling's formula [2 5) to approximate the gamma function. 



exp 



r(5Z/2) 

Next we expand the exponent around the saddle point 



/(Z)^(Z-l)lne- 



5Z 1\ , 5Z 5Z 1 , , 
^-:^ ln^ + ^- :^ ln(27r) 
2 2/2 2 2^^ 



f{Z* 



-f"{Z*){Z-Z*f 



where Z* is the saddle point value of Z, which satisfies the saddle point equation. 

/'(Z*)=lne--ln— + — = 
J \ ) ^2 2 2Z* 

The distribution function (|43p can be approximated as a Gaussian. 

exp [nZ*){Z - Z*)y2] 



Pc^iZ) 



/OG 
dZexp [f"{Z*){Z-Z*f/2\ 
-oo 



bZ* 



The average value of Z is calculated to be 

{Z)cq_ = Zq 



hZ* 



1 



4(Z* 



{Z -Z* 



(102) 



(103) 



(104) 



(105) 



(106) 



Therefore we can replace Z* in eq (jlOSp by Zq. Finally we have eq (|44l) as the approximate form 
for the subchain number distribution function. 



C.2 Saddle Point Approximation for Pcq('^) 

We consider to obtain a simple approximate expression for eq (jSOp . For sufficiently large Zq, we 
can use the following approximate form. 
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and Z** is given via the following saddle point equation. 

g'{Zn = -ln(l--) +ln£--ln — 
^^^2V nJ ^ 2 2 



1 
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(107) 



From eqs (|104p and (|106p . we have the following approximate relation between Z** and Zq 

n ' 
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N 
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(108) 



(109) 



(110) 
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Then we can write 



Z=2 



2. TTTTTT^^TT:^^ — (l--)exp 



r(5(z-i)/2) 
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2 V iV 



(111) 



and the segment number distribution function can be expressed approximately as 



3/2 



3 /2 

C!C n ' exp 



5n 
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1 - ^ ) exp 



5^0 / _ n_\ 
2 V iV/ 



(112) 



By normalizing eq (|112p . finally we have eq ((5T|) as the approximate form for Pcq(n). 

D Strong Repulsion Limit 

In the main text and Appendix [Aj we calculated some equilibrium distribution functions for re- 
pulsive and non-interacting slip-link models, respectively. For comparison, here we consider the 
case where the repulsive interaction between slip-links is quite strong. We call such a limit as the 
strong repulsion limit. As we will show, slip-links form a Wigner crystal like ordered structure in 
this limit. 

We assume the following effective potential for slip-links, instead of eq (PU)) . 



C^siip-iink("-) = -afcsTlnn 



(113) 



Here, a is the parameter which represents to the strength of the repulsive interaction, (a — and 
a = 3/2 corresponds to the non- interacting and repulsive slip-link models, respectively.) At the 
strong repulsion limit, we set a — > oo. 

The equilibrium probability distribution function and the grand partition function are expressed 
as follows. 



P,q{{Rk},{Nk},Z) = - 
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(114) 
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3{Rk — Rk-i) 
2iVfc62 



(115) 

Since the repulsive interaction energy becomes very large at the strong repulsion limit, we can 
reasonably utilize the saddle point approximation for {Nk}- For a ^ 1, the distribution of {Nk} 
is expected to be sufficiently sharp and thus the fluctuations around the saddle point is negligible. 
Then the grand partition function can be approximated as 



^^z+D^z-i / d{Rk} ( ^ 
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(116) 



For simplicity, we assume that the average number of subchains is sufficiently large (Zq 3> 1) 
and utilize the saddle point approximation for Z . The result is 



( ^ 



\2ttZo 
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exp 



2Z, 
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(117) 



5{Z — Zq) (a — > oo) 
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Thus we find that the subchain number distribution simply becomes the delta function. This means 
that there is essentially no subchain number fluctuation in the strong repulsion limit. Similarly, 
the segment number distribution function is obtained as 

Peq(n)«5(n-iVo) (118) 

Eq (jllSp means that the chemical distance between two neighboring slip-links is constant (slip-links 
are placed equidistantly on a chain) . Eqs (|117p and (|118p are naturally derived from the dependence 
of the variances of Z and n on a. From the results in Appendix|3 the main text, and this appendix, 
the variances are roug hly estimated as ((Z-Zo)^)oq~ Zo/(l + a) and ((n-7Vo)^)cq ^Nll{\^a). 
At the strong repulsion limit (a — >■ 0), both of them approach to zero and the distribution functions 
reduce to delta functions (as eqs (|117p and (jllSp ). We can interpret such a state as a sort of Wigner 
crystal. This situation seems to be qualitatively different from most of slip-links models such as the 
PCN model. It would be rather similar to the simple tube model [5] in which each tube segment 
contains a constant number of segments. (But of course, the slip-link model is not equivalent to 
the tube model.) Anyway, the distribution functions for Z and n at the strong repulsion limit 
are much sharper than distribution functions for the non-interacting or repulsive slip-link models. 
Thus we consider that these distribution functions become sharper as the repulsive interaction 
between slip-links increases. 

Because the segment number n is fixed to A^Oi the bond vector distribution function is nothing 
but a Gaussian. The bond vector and bond length distribution functions become 

At the limit of small Q, we have Pcq(Q) cx as the asymptotic form. This asymptotic behavior is 
similar to one of the repulsive slip-link model (eq (|59|) ). Eq (|120p has a maximum at QI\fN^ — 
■y/2/3 w 0.816, which is larger than the value for the repulsive slip-link model. Thus we expect that 
the value of Q/ \/NqE^ which gives the maximum increases as the repulsive interaction increases. 
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Figure Captions 



Figure [TJ The slip-linked subchain number distribution functions with various values of Zq. Solid 
and dashed curves show the exact and approximate distribution functions for the single chain 
repulsive slip-link model, respectively. The dotted curves show the Poisson distribution which cor- 
responds to the distribution of the single chain non-interacting slip-link model. Zq = 5, 10, 20, 40, 
and 80 from left to right. The approximate distributions almost coincide to the exact distributions. 

Figure[2j Dependence of the variance of slip-linked subchain number {{Z — 2^o)^)cq on the average 
{Z)eq = Zq. The solid and dashed curves show the exact and approximate variances for the 
repulsive slip-link model, and the dotted curve shows the variance for the non-interacting slip-link 
model {{{Z — 2'o)^)oq = Zq — 1). The approximate curve is very close to the exact curve. 

Figure [31 The segment number distribution functions for sufficiently large Zq. Solid and dot- 
ted curves show the distribution functions for the repulsive and non-interacting slip-link models, 
respectively. 

Figure m The bond length distribution functions for sufficiently large Zq. Solid and dotted curves 
show the distribution functions for the repulsive and non-interacting slip-link models, respectively. 

Figure [5j The slip-linked subchain number distribution functions by single chain models and PCN 
simulations. Curves are theoretical predictions by single chain models and symbols are simulation 
data. Zq = 5, 10, 20, and 40 from left to right. 

Figure m The segment number distribution functions calculated by single chain models and PCN 
simulations (for Zq = 40). Curves and symbols represent the results by single chain models and 
PCN simulations, respectively. 

Figure [T] The bond length distribution functions calculated by single chain models and PCN 
simulations (for Zq — 40). Curves and symbols represent the results by single chain models and 
PCN simulations, respectively. 
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